%% Generate Random Rotation
load QQ_full % Load pre-saved seeds used in the paper
% QQ =[];
% NQ= 1.5*10^(6); %Number of Random Rotation
% parfor iii = 1:NQ
%     %generate a random rotation
%     v = randn(VAR.neqs,VAR.neqs); %First Generate a totally random 3x3 matrix v
%     [q, r]=qr(v,0); %Obtain the orthonormal matrix of v, call it q
%     QQ(:,:,iii)=q*diag(sign(diag(r)));    % this makes the diagonal of r positive
% end
% save QQ_full QQ
yyyy=repmat((1961:2019)',1,12);
yyyy=reshape(yyyy' ,size(yyyy,1)*12,1);
mm  =repmat((1:12)',size(yyyy,1)/12,1);
dates_lag=yyyy*100+mm;
%datem=datem(7:end); 
 dates_lag=dates_lag(VAR.nlag+1:end); %
clear yyyy mm
%% Compute e
 %Ebar = [3.8720 3.7772]; % see page 21 "The 75th-percentile value of eFt  in 1987:10 is 4.16, and is 4.57 in 2008:09"     
 %Ebar = [3.5894 4.4519]; %  with 4 variables (UNC, IP , employment, STATE = Unemployment)
  Ebar = [3.5462 3.4036]; %  with 4 variables (UNC, IPG, employment, STATE = Unemployment)
% % Generate Ebar
% Ebar                         = gen_k(VAR.residuals, QQ,                            dates_lag,     pick);

load ExtInstr
%plot([diff(VAR.Y(:,2)) ExtInstr(VAR.nlag+2:end,1)])
%corr([diff(VAR.Y(:,2)) ExtInstr(VAR.nlag+2:end,1)])
% Generate B
B_out                          = gen_B(VAR.residuals,   QQ,ExtInstr(VAR.nlag+1:end,:),dates_lag,Ebar,pick);
% Gen ehat
[B_allsol,max_Ck_indx]     = gen_lo_hi(VAR.residuals,B_out,ExtInstr(VAR.nlag+1:end,:),               pick);
Bfinal = B_allsol(:,:,max_Ck_indx);
save Bfinal Bfinal
 %ehat_maxG = VAR.residuals*inv(Bfinal)';
  ehat_maxG =                  (Bfinal\VAR.residuals')';
%#Check
% ehat_maxG(index_87       ,pick)
% ehat_maxG(index_lehman   ,pick)
% ehat_maxG(index_79       ,pick)
% min(ehat_maxG(index_crsis,pick))
% corr(  Vs(VAR.nlag+1:end), ehat_maxG(:,pick))
% corr(Gold(VAR.nlag+1:end), ehat_maxG(:,pick))